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ABSTRACT 


Electrical  distribution  systems  are  the  last  stage  of  electricity  delivery,  bridging  the  gap 
between  transmission  networks  and  the  customer.  Compared  to  transmission  systems, 
there  has  been  comparatively  little  analysis  performed  on  distribution  infrastructure.  We 
first  establish  a  minimum  load  shed  optimal  switching  model.  The  model  seeks  to  ensure 
all  electricity  demand  is  met  for  a  given  distribution  system  by  minimizing  the  amount  of 
load  not  met.  This  is  accomplished  by  closing  and  opening  switches  within  the  network, 
subject  to  various  electrical  physics  constraints.  We  initially  test  our  model  against  a  simple 
toy  network  to  verify  proper  operation  before  expanding  to  increasingly  complex  networks. 
We  culminate  with  testing  against  a  standardized  test  circuit. 

The  model  runs  in  an  open-source  programming  language  on  a  consumer-grade 
lap-top,  reading  in  data  from  common  formats  with  solvers  readily  familiar  to  the 
optimization  community. 
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Executive  Summary 


Of  the  16  critical  infrastructure  sectors,  electric  power  is  identified  as  being  significantly 
important  because  it  enables  the  operation  of  the  other  sectors.  Within  the  realm  of  electric 
power,  electrical  distribution  systems  consist  of  the  substations,  transformers,  and  lines  that 
bridge  the  gap  between  transmission  systems,  which  originate  at  power  plants  such  as  coal 
or  nuclear,  and  the  end-user.  Distribution  systems  are  most  readily  characterized  as  the 
power  poles  and  lines  running  between  homes,  commercial  businesses,  and  neighboring 
municipalities. 

Working  in  concert  with  our  associates  at  Los  Alamos  National  Labs  (LANL)  and  under 
the  Joint  Warfare  Analysis  Center  (JWAC)  -  NPS  research  program,  we  build  a  linearized 
optimization  model  in  the  Python  programming  language,  implementing  the  mathematical 
program  in  the  Pyomo  environment.  The  original  formulation  was  provided  by  our  LANL 
associates  and  expanded  upon  during  model  development.  At  its  core,  the  model  delivers 
real  and  reactive  power  from  an  originating  substation(s)  to  loads  throughout  the  network, 
adhering  to  various  electrical  engineering  physics  constraints. 

These  constraints  include  flow  balance,  balance  between  the  power  transmission  for  phases 
within  a  single  line,  adhering  to  the  thermal  limits  of  lines,  line  losses  due  to  resistance  and 
reactance,  maintaining  nominal  system  operating  voltage,  and  various  topology  constraints. 
The  entire  model  is  driven  by  an  objective  function  seeking  to  minimize  the  load  shed 
throughout  the  network. 

Testing  and  analysis  is  initially  performed  on  a  test  circuit  provided  by  our  associates  at 
LANL.  The  circuit  is  simple  enough  that  basic  operation  can  be  checked  by  hand  calculation 
to  verify  model  behavior.  We  then  specifically  examine  behavior  of  the  topology  constraints 
against  another  test  circuit  we  designed.  We  observe  the  model  performing  proper  operation 
of  switches  within  the  network  to  preclude  cycle  formation  while  still  minimizing  load  shed. 

Our  testing  culminates  with  application  of  an  Institute  of  Electrical  and  Electronics  Engi¬ 
neers  (IEEE)  standardized  test  circuit.  The  test  circuit  possesses  more  robust  and  dynamic 
electrical  characteristics,  specifically  imbalanced  loading  and  phase  variations.  Zero  load 
shed  is  expected  with  the  test  circuit  and  we  observe  none  with  our  model.  The  downstream 


xv 


goal  of  our  work  here  is  to  perform  analysis  on  real-world  distribution  systems  when  they 
are  subject  to  various  disruptions.  As  an  intermediate  step  toward  that  goal  we  duplicate 
the  IEEE  network  into  “a”  and  “b”  sides,  connecting  them  with  a  single  line  containing  a 
switch.  We  proceed  to  remove  the  substation  from  the  “b”  side  and  observe  the  resulting 
behavior  as  the  remaining  “a”  substation  attempts  to  power  the  entire  system. 

We  conclude  with  recommendations  for  continued  work  and  analysis.  Specifically,  our 
model  can  be  folded  into  attacker-defender  analyses  to  inform  researchers  on  vulnerabilities 
within  a  given  distribution  system.  Furthermore,  testing  of  real-world  systems  by  triggering 
systematic  failures  of  given  components  can  yield  insights  to  distribution  system  behavior 
when  under  duress.  Our  work  is  unique  due  to  its  focus  on  distribution  systems,  its  adherence 
to  electrical  engineering  mathematics,  and  results  with  guaranteed  optimality. 
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CHAPTER  1: 
Introduction 


The  U.S.  government  defines  critical  infrastructure  as  “systems  and  assets,  whether  physical 
or  virtual,  so  vital  to  the  United  States  that  the  incapacity  or  destruction  of  such  systems  and 
assets  would  have  a  debilitating  impact  on  security,  national  economic  security,  national 
public  health  or  safety,  or  any  combination  of  those  matters”  (White  House,  2013).  Critical 
infrastructure  systems  are  those  which  deliver  the  energy  to  heat  and  cool  our  homes,  the 
water  we  drink,  the  communication  systems  we  enjoy,  and  the  transportation  systems  that  are 
integral  to  the  American  way  of  life.  Furthermore,  these  systems  serve  as  the  foundation  for 
the  economy,  security,  and  health  of  our  nation  (Department  of  Homeland  Security,  2017). 
Due  to  the  significance  of  these  systems,  there  is  considerable  interest  in  understanding, 
modeling,  and  predicting  of  infrastructure  system  response  to  disruption.  It  is  of  paramount 
interest  at  the  presidential  level  that  we  “maintain  secure,  functioning,  and  resilient  critical 
infrastructure”  (White  House,  2013).  Electric  power,  one  of  these  critical  infrastructure 
systems,  is  the  focus  of  this  research. 

1.1  The  Importance  of  Electric  Power 

As  defined  by  the  Department  of  Homeland  Security  (DHS),  there  are  16  critical  infras¬ 
tructure  sectors  (Department  of  Homeland  Security,  2017).  Of  these,  electric  power  is 
particularly  important  because  it  provides  an  “enabling  function”  for  the  others  in  funda¬ 
mental  ways  (White  House,  2013). 

The  overall  structure  of  electrical  power  infrastructure  can  be  loosely  categorized  into 
three  main  parts:  generation,  transmission,  and  distribution.  Generation  systems,  such  as 
solar,  hydro,  wind,  nuclear,  or  gas  turbine  power  plants,  produce  electricity.  Transmission 
systems  move  electricity  in  bulk  from  the  originating  power  plant  to  substations  and  are 
most  commonly  recognized  as  the  massive  towers  carrying  heavy  lines  across  rural  stretches 
of  highway.  The  distribution  system  exists  between  the  substation  and  end  customer,  most 
readily  characterized  as  the  telephone  poles  abundant  in  any  residential  and  commercial 
area. 
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Disruptions  to  the  electric  power  system  can  be  caused  by  a  variety  of  sources.  Foremost 
among  these  are  natural  disasters  such  as  hurricanes  or  earthquakes.  For  example,  Hurricane 
Sandy  in  2012  caused  power  losses  to  over  8,000,000  homes  across  17  states  (Webley,  2012). 
Other  potential  causes  of  electric  power  disruption  are  operator  error  (e.g.,  the  Three  Mile 
Island  incident;  see  Nuclear  Regulatory  Commission,  2013)  or  the  actions  of  an  adversary 
(e.g.,  the  coordinated  attacks  on  the  Metcalf  substation  in  2013  near  San  Jose,  CA;  see 
Smith,  2014).  Consequences  for  these  disruptions  vary  depending  on  the  scope  and  scale  of 
the  damage  as  well  as  the  types  of  electricity  infrastructure  damaged.  A  car  hitting  a  power 
pole  may  interrupt  power  to  a  single  home,  but  a  worst-case  solar  event  electromagnetic 
pulse  (EMP)  would  take  an  estimated  four  to  ten  years  to  repair  and  leave  most  of  the  United 
States  without  power  (Federal  Energy  Regulatory  Commission,  2010). 

Regardless  of  the  cause,  the  consequences  of  a  disruption  have  significant  negative  impli¬ 
cations  for  our  welfare,  both  in  terms  of  loss  of  life  as  well  as  negative  economic  impacts 
(Anderson  and  Bell,  2012;  Rose  et  al.,  1997).  Generation  facilities  are  frequently  taken 
online  and  offline  to  balance  power  supply  and  load  demand,  but  unexpected  failure  of  a 
generation  facility  can  cause  problems  (Knaus,  2017).  Disruptions  to  transmission  systems, 
either  from  the  loss  of  a  high-voltage  line  or  a  substation,  can  impede  the  ability  to  move 
bulk  power  and  cause  long-distance  power  outages  (Cava,  2017).  A  disruption  to  a  distri¬ 
bution  system  causes  downstream  customers  to  potentially  go  without  electricity,  but  these 
effects  are  often  local  in  nature. 

There  is  a  tactical  need  to  understand  the  strengths  and  weaknesses  of  our  own  electric  grid 
as  well  as  those  of  potential  adversaries.  Domestically,  identifying  weaknesses  in  our  grid 
informs  our  defense  efforts  so  we  can  minimize  the  impact  of  adversarial  action  against 
us,  because  as  reporters,  politicians,  and  utility  operators  know,  disruptions  are  devastating 
(Williams  and  Bennett,  2016).  From  an  offensive  perspective,  we  desire  to  know  the  effects 
of  targeting  an  enemy’s  infrastructure  to  achieve  our  tactical  objective,  but  we  desire  to 
restrict  the  consequence  to  only  what’s  necessary  to  achieve  victory.  Militarily,  the  second 
order  effects  of  taking  out  an  entire  transmission  system  in  an  enemy  territory  might  prove 
too  costly.  Research  has  shown  that  in  large  scale  power  outages,  ordinary  civilians  suffer 
and  death  rates  increase  (Anderson  and  Bell,  2012).  If  we  seek  to  confine  power  disruption 
to  as  small  an  area  as  possible,  attacking  the  transmission  system  is  off  the  table.  Hart  et  al. 
(2014)  provide  a  discussion  of  the  tensions  and  tradeoffs  for  targeting  infrastructure  in  a 
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military  context. 


Prior  research  and  analysis  efforts  at  the  Naval  Postgraduate  School  (NPS)  have  focused 
primarily  on  the  behavior  and  response  of  transmission  systems.  There  are  multiple  re¬ 
search  papers  and  theses  on  the  subject  of  transmission  system  modeling  (e.g.,  Rose,  2007; 
Salmeron  et  al.,  2012).  Comparatively,  little  effort  has  been  applied  to  distribution  systems. 

The  focus  of  this  thesis  is  understanding  how  a  distribution  system  will  respond  to  the  loss 
of  components  (e.g.,  substations),  and  in  particular  the  presence  and  size  of  outages.  In 
response  to  component  failure,  a  distribution  system  will  adjust  its  connectivity  by  opening 
and  closing  switches  to  continue  powering  as  much  of  the  grid  as  possible.  We  capture  this 
behavior  in  a  prescriptive  decision  Alternating  Current  Power  Flow  (ACPF)  model  which 
changes  flow  by  changing  connectivity  via  switch  position  adjustments  to  minimize  the 
amount  of  load  shed. 

We  implement  the  ACPF  model  in  the  Python  programming  language  (Python  Software 
Foundation,  2017),  primarily  relying  on  the  Pyomo  optimization  package  (ffart  et  al.,  2012, 
201 1).  We  verify  model  operation  by  application  to  a  toy  model  and  then  expand  to  a  real 
data  set  from  a  large  Northeastern  utility  provider. 

1.2  Thesis  Structure 

Chapter  2  provides  cursory  overviews  of  the  academic  disciplines  underpinning  this  work. 
Additionally,  we  review  the  relevant  literature  and  research  pertaining  to  our  analysis.  In 
Chapter  3  we  introduce  the  ACPF  model  and  demonstrate  its  application  with  a  small  series 
of  toy  data  sets.  Chapter  4  describes  our  application  of  the  model  to  a  standardized  test 
circuit.  We  conclude  in  Chapter  5  by  summarizing  our  findings  along  with  recommendations 
for  future  research. 
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CHAPTER  2: 
Background 


This  chapter  provides  background  information  on  the  disciplines  involved  in  this  thesis  as 
well  as  pertinent  previous  research. 

2.1  Technical  Background  Overview 

Our  work  encapsulates  two  major  academic  domains:  mathematical  optimization  and 
Electrical  Engineering  (EE).  While  it  is  beyond  the  scope  of  this  thesis  to  fully  inform 
the  reader  on  both  these  subjects,  we  seek  to  provide  a  brief  overview  to  aid  reading  our 
work. 

2.1.1  Mathematical  Programming 

Mathematical  programming  involves  theory  and  methods  for  solving  for  extreme  points, 
either  maxima  or  a  minima,  of  functions  on  sets  defined  by  linear  and/or  non-linear  con¬ 
straints.  The  primary  objective  function  for  which  a  maxima  or  a  minima  is  being  found, 
commonly  involves  a  maximum  profit  or  lowest  cost.  The  constraints  are  linear  equality 
and  inequalities  representing  limitations  imposed  on  the  system,  such  as  how  much  money 
can  be  spent  or  the  capacity  limits  of  a  path.  For  a  thorough  introduction  to  optimization 
and  mathematical  programming,  see  the  textbook  Optimization  in  Operations  Research 
(Rardin,  2016). 

2.1.2  Electrical  Engineering 

For  a  deeper  understanding  of  EE  as  it  applies  to  electric  power  systems  we  recommend  the 
work  Elements  of  Power  Systems  (Sadhu  and  Das,  2015). 

Electrical  transmission  and  distribution  systems  are  frequently  3-phase  Alternating  Current 
(AC)  systems.  The  three  phases  are  conventionally  a,  b,  or  c  and  a  given  load  or  line  in  the 
network  can  be  any  combination  of  the  three  phases.  The  supply  and  demand  quantities 
in  the  network  are  measured  in  watts,  with  transmission  systems  generally  operating  in 
Megawatts  (Mw)  while  distribution  systems  frequently  are  measured  in  Kilowatts  (Kw). 
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Generated  electricity  consists  of  real  and  reactive  components.  Real  power  produces  kinetic 
work  (e.g.,  motors  spinning),  while  reactive  power  is  the  electromagnetic  hindrance  of 
electrical  processes.  Correspondingly,  loads  in  the  network  consume  real  and/or  reactive 
power.  Individual  power  lines  possess  both  resistance  and  reactance,  qualities  that  hinder 
real  and  reactive  power  flow,  respectively. 


2.2  Previous  Research 

2.2.1  Attacker-Defender  and  Defender- Attacker-Defender  Analysis 

Attacker-Defender  (AD)  analysis  is  a  two-part  optimization  problem.  The  attacker  seeks  to 
maximize  the  minimum  value  the  defender  can  obtain  from  the  system  post  attack.  There 
are  three  key  assumptions  in  AD  modeling,  “(1)  the  attacker’s  and  defender’s  actions  are 
sequential,  (2)  the  attacker  has  a  perfect  model  of  how  the  defender  will  (or  should)  optimally 
operate  the  system,  even  after  an  attack,  and  (3)  the  attacker  will  manipulate  that  system  to 
his  best  advantage”  (Brown  et  al.,  2006). 

Defender  is  a  bit  of  a  misnomer  in  an  AD  model.  Nothing  is  being  hardened  or  made 
invulnerable  to  attack;  in  reality  the  defender  is  simply  a  system  user  who  best  operates  post¬ 
attack.  However,  it  is  possible  to  expand  a  two-level  AD  model  into  a  three-level  Defender- 
Attacker-Defender  (DAD)  model,  where  the  system  is  initially  hardened  to  preclude  or 
mitigate  the  most  damaging  attack  (Brown  et  al.,  2006).  Afterward  and  with  knowledge  of 
the  defensive  action,  the  attacker  performs  a  most  damaging  attack,  and  lastly  the  system 
operator  works  to  obtain  best  system  performance. 

Salmeron  et  al.  (2004)  applies  attacker-defender  optimization  supported  by  a  heuristic  to 
an  electric  power  model.  They  account  for  a  great  deal  of  electrical  physics  considerations, 
abstaining  from  the  simplifications  of  a  pure  flow  model.  Their  results  demonstrate  the 
ability  to  identify  the  key  components  in  a  network  that,  if  attacked,  would  be  most  crippling 
to  operation  and  simultaneously  represent  the  greatest  need  for  hardening  and  resiliency 
efforts.  Analysis  is  performed  on  Institute  of  Electrical  and  Electronics  Engineers  (IEEE) 
Reliability  Test  System  (RTS)  96. 

Salmeron  et  al.  (2005)  introduce  a  Graphical  User  Interface  (GUI)  for  electric  power 
infrastructure  analysis,  named  Vulnerability  of  Electric  Grids  Analyzer  (VEGA),  which 
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is  supported  by  an  optimal  power  flow  linear  program  and  an  interdiction  model  based 
on  heuristic  algorithms.  The  heuristic  algorithm  is  used  to  “find  the  combination  of 
facilities  (generating  units,  buses,  transmission  lines,  transformers  and/or  substations)  whose 
destruction  would  cause  the  largest  or  nearly  largest  disruption  in  the  system  at  a  given 
point  in  time,  which  in  turn  suggests  the  system  components  most  in  need  of  hardening” 
(Salmeron  et  al.,  2005).  Tests  and  analysis  are  performed  on  two  real  world  transmission 
systems.  Additionally,  a  time-series  component  is  incorporated  with  acknowledgement  of 
peak  loading  time-frames  being  the  most  devastating  windows  of  attack. 

Salmeron  and  Wood  (2009)  extend  prior  work  to  formulate  and  solve  DAD  problems  for 
electrical  transmission  grids,  identifying  critical  components  sensitive  to  terrorist  attack. 
Applying  a  generalized  version  of  Benders  decomposition,  Salmeron  et  al.  (2009)  identify  a 
set  of  components  whose  interdiction  causes  maximum  economic  loss  to  customers.  Their 
test  problems  “describe  a  regional  power  grid  in  the  United  States  with  approximately  5000 
buses,  6000  lines,  and  500  generators”  (Salmeron  et  al.,  2009). 

Salmeron  et  al.  (2011)  perform  DAD  analysis  on  Vandenberg  Air  Force  base  and  its  sur¬ 
rounding  area.  Their  efforts  identify  the  vulnerable  components  whose  destruction  would 
be  maximally  disruptive  to  system  operation  and  presume  adversaries  are  aware  of  those 
defensive  preparations,  in  turn  making  maximally  destructive  attacks.  Lastly,  the  system 
is  optimally  operated  post-attack.  This  focus  on  a  Department  of  Defense  (DOD)  entity  is 
unique  because  the  surrounding  area  benefits  from  excess  generation  capacity. 

Further  modeling  efforts  focus  on  the  Guam  Power  Authority  (GPA),  identifying  critical 
sets  of  components  whose  destruction  would  most  impede  GPA’s  ability  to  provide  power 
(Salmeron  et  al.,  2012).  Subsequently,  these  sets  of  components  benefit  the  most  from 
hardening  efforts. 

2.2.2  Previous  NPS  Theses 

Analysis  on  electric  power  systems  using  AD  and  DAD  models  has  been  driven  and  sup¬ 
ported  by  many  student  theses. 

Stathakos  (2003)  performs  research  contributing  to  the  development  of  VEGA.  Alvarez 
(2004)  demonstrates  empirical  speed  ups  to  Benders  decomposition  on  electric  grid  inter- 
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diction  problems,  tested  against  IEEE  RTS  networks.  Optimality  gap  reductions  occurred 
alongside  speed  improvements. 

Carnal  (2005)  further  improves  on  VEGA  by  inclusion  of  the  Xpres-MP  software  package, 
reducing  run-time  by  65%  -  85%. 

Schneider  (2005)  performs  DAD  analysis  on  the  eastern  interconnection  transmission  grid 
and  determine  hardening  levels  necessary  to  reduce  risk  below  acceptable  levels.  Jonart 
(2008)  does  follow-on  work,  making  algorithmic  improvements,  yielding  results  within  1% 
of  optimality. 

Ang  (2006)  performs  a  Mixed-integer  program  (MIP)  optimization  of  electric  grid  recovery 
efforts  after  natural  disaster  or  terrorist  attack.  Constraints  include  availability  of  repair 
resources,  penalties  for  unserved  electricity  demand,  and  repair  time  planning  horizons. 
Testing  is  performed  on  an  IEEE  300  bus  transmission  network. 

Rose  (2007)  develops  a  mathematical  model  to  identify  critical  components  to  defend  in  an 
electric  grid,  given  limited  resources.  Using  a  pure-flow  model  for  the  electric  grid.  Rose 
reduces  solution  times  from  previous  efforts  while  simultaneously  finding  solutions  that  are 
40%  better. 

2.2.3  Damage  Anticipation  Modeling  on  Electrical  Systems 

Guikema  et  al.  (2010)  use  regression  and  data  mining  techniques  to  estimate  required 
repair  time  and  replacement  utility  poles  and  transformers  in  the  wake  of  electric  power 
infrastructure  damage  resulting  from  natural  disaster,  notably  hurricanes.  Their  efforts 
are  predictive  in  nature  and  seek  to  best  prepare  responders  for  repairs  after  a  disruption. 
Ouyang  and  Duenas-Osorio  (2014)  build  a  probabilistic  model  to  determine  economic, 
social,  technical,  and  organizational  resilience  to  hurricane  damage. 

Maliszewski  and  Perrings  (2012)  use  regression  modeling  to  build  spatial  outage  visual¬ 
izations  and  identify  geographic  factors  with  strong  correlation  to  minimal  outage  risks  in 
Phoenix’s  electrical  distribution  system.  They  find  a  strong  correlation  of  neighborhoods 
near  high  priority  customers  such  as  hospitals  experiencing  shortened  outages  because  of 
“the  quality  of  infrastructure  provided  and  the  priority  given  to  restoring  supply  to  those 
users”  (Maliszewski  and  Perrings,  2012).  Proximity  to  vegetation  was  statistically  signifi- 
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cant  as  well  while  other  services,  such  as  police  and  fire  stations  or  schools,  all  which  were 
insignificant  or  negatively  correlated. 

A  variety  of  natural  disasters,  notably  the  2010  earthquakes  in  Haiti  and  Chile,  the  2011 
earthquake  and  tsunami  in  Japan,  as  well  as  multiple  hurricanes,  have  served  as  the  sources 
of  data  for  various  resiliency  modeling  efforts.  MacKenzie  and  Barker  (2012)  build  a  re¬ 
gression  model  and  mixed  effects  model  to  anticipate  rates  of  recovery  post-disruption  and 
anticipate  a  system’s  resiliency.  Using  data  from  Hurricane  Ike,  Wei  et  al.  (2013)  develop  a 
“spatial-temporal  non- stationary  random  process  to  model  large-scale  disruptions  of  power 
distribution  induced  by  severe  weather.”  Their  stochastic  approach  yields  resiliency  esti¬ 
mates  for  components  and  identifies  dynamic  neighborhoods  within  the  network.  Kwasinski 
(2016)  performs  data  analysis  and  regression  modeling  to  inform  future  damage  anticipation 
efforts  and  design  considerations  for  distribution  systems. 

2.2.4  Optimal  Switching  and  RTS  usage 

Delgadillo  et  al.  (2010)  construct  an  optimization  approach  to  defend  a  grid  against  terrorist 
attack  by  calculating  the  best  operator  response  (e.g.,  switching)  to  disruption.  Their 
work  is  a  two-part  optimization,  first  identifying  the  most  damaging  terrorist  attack, 
measured  in  Mw  not  delivered,  and  then  identify  the  optimal  switch  realignment  in 
response  to  that  attack. 

Additionally,  de  Assis  et  al.  (2015)  consider  optimization  of  switch  allocation  for  distribution 
systems.  They  use  a  distribution  system  and  identify  the  ideal  locations  for  switches  under 
various  system  design  constraints.  Specifically,  their  research  minimizes  the  cost  of  switch 
allocation,  a  system  design  consideration.  They  also  minimize  energy  not  supplied  while 
meeting  various  regulatory  agency  reliability  requirements  and  system  capacity  necessities. 
Their  solutions  are  determined  by  “memetic  algorithms,  supported  by  well-established 
genetic  algorithm  concepts  and  local  search  optimization  procedures.” 

The  IEEE  118-bus  and  the  RTS  96  test  system  are  studied  for  optimal  transmission  system 
switching  with  a  mixed  integer  program  (Hedman  et  al.,  2009).  Their  model  demonstrates 
the  test  networks  can  be  operated  to  meet  N-l  secure  standards  while  simultaneously 
cutting  costs.  Another  instance  of  RTS  usage  is  the  work  of  Alguacil  et  al.  (2003)  in  their 
transmission  expansion  planning  optimization  problem. 
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CHAPTER  3: 
Model 


We  develop  an  alternating  current  (AC)  model  for  electrical  distribution  system  that  in¬ 
tegrates  switching  and  adherence  to  full  3-phase  AC  mathematics.  Instead  of  suggesting 
likelihoods  of  what  could  happen  in  a  disaster  or  terrorist  scenario,  we  model  what  will 
happen  for  any  disruption  or  series  of  disruptions. 

3.1  Model  Formulation 

Our  Alternating  Current  Power  Flow  (ACPF)  model  is  built  in  collaboration  with  our 
partners  at  Los  Alamos  National  Laboratory  (LANL),  notably  Harsha  Nagarajan  and  Arthur 
Barnes,  with  significant  sourcing  from  their  previous  work  (e.g.,  Barnes  et  al.,  2017; 
Nagarajan  et  al.,  2016). 

The  entirety  of  our  ACPF  model  is  presented  in  NPS  format  for  a  given  Power  Grid  (PG). 
Discussion  of  the  objective  function  and  constraints  follows  in  section  3.1.2. 

3.1.1  Mathematical  Formulation 

Sets  and  Indices 


i  g  N 
ScNxN 
NT  c  <5 
SQ  NT 
3K  c  NxN 

P  -  {a,  b,  c} 


set  of  nodes  (buses)  in  PG  (alias  j) 

set  of  undirected  edges  (lines  and  transformers)  in  PG 

set  of  edges  without  transformers  in  PG 

set  of  edges  with  switches  in  PG 

set  of  directed  arcs  (lines  and  transformers)  in  PG 

(l  j)  e  £  (6  j)  6  £  A  (j,  i)  G  £ 

set  of  phases  in  PG 
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k  e  Pi  c  <P  set  of  phases  on  bus  i  which  contains  a  load  or  a  generator 


Pij  c  P  set  of  phases  on  line  (ij)  in  PG 
C  set  of  cycles  in  PG 
N(C)  set  of  nodes  in  cycle  C 
<5(C)  set  of  edges  in  cycle  C 

Parameters  [units] 

Api  k  cost  for  shedding  real  power  at  bus  i  on  phase  k  [$/watt] 

Aqi  k  cost  for  shedding  reactive  power  at  bus  i  on  phase  k  [$/watt] 

Tjj£  thermal  limit  on  edge  (i.j)  for  phase  k  [Amps] 
vj , ,  vu,  voltage  lower  and  upper  bounds  at  bus  i  for  phase  k  [Kv] 
i  imaginary  number  constant 

dpj  k  +  i  dcp  k  AC  power  demand  at  bus  i  on  phase  k  [Kw  and  Kvar] 

gpuik  +  i  gq1' k  existing  AC  power  generation  capacity  at  bus  i  on  phase  k  [Kw  and  Kvar] 

Rjj'k  +  i  Xjj,k  impedance  of  edge  ( i,j )  on  phase  k  in  PG  [Q] 

Zij  impedance  of  edge  (IJ).  defined  as: 

Z,J  =  z,je^-,  Zu  =  Zije-W  where,  Z,j  =  R,j  +iXlj.  Z,,  =  R„+i XtJ  [01 
Mu  Mij  =  (1.05V”OW)2  -  (0.95 P""'")2  V/./  6  6  [Volts2] 

Pg  allowed  variation  in  flow  between  phases  on  edge  (ij)  [0.15  for  transformer 
edges,  1.0  for  all  other  edges] 

yn°m  nomjna]  voltage  the  system  is  operating  at  [Kv] 
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Decision  Variables  [units,  when  appropriate] 


vik  voltage  at  bus  i  on  phase  k  [Kv] 

Xij  determines  if  directed  arc  (i,j)  carries  power  flow  in  PG  [binary] 
Pij  k  real  AC  power  flow  on  arc  (i,  j )  [Kw] 
qjhk  reactive  AC  power  flow  on  arc  (i,j)  [Kvar] 

Ipi'U  +  i  lqi,k  AC  power  load  that  is  shed  at  bus  i  on  phase  k 


gpi^k  +  i  gqi,k  existing  AC  power  generation  capacity  at  bus  i  on  phase  k 


Formulation 

min  I  (Api,klPi,k  +  Mi,klquk )  (1) 

i£N,k£Pj 

S.t. 


Power  flow  constraints 


9Pi,k  -  dpi,k  +  lpi,k  =  J]  Pij,k  -  J]  PiUk 

i jeS  ji&£ 

gqik  —  dqi  k  +  lqi,k  —  ^  ^  qij,k  ~  ^  ]  q ji.k 

ij£&  ji£& 


(l-Av) 


YjkePij  Pij,k 


<  Pij,k  ^  (1  +  fiij) 


JukePi,  Pij,k 


—  qij,k  —  ( 1  +  Ay) 


TjkePtj  qij,k 


Vi  €  ft,  k  £  Pi 


Vi  £  ft,  k  £  Pj 


Vi  j  £  ft,  k  e  Pij 


Vi  j  £  ft,  k  £  Pij 

Vi  j  £  ft,  k  £  Pi 
Vi  j  £  ft,  k  £  Pi 


a  -  Ai)  s  S(U  Ay)  v,7  ao  Ay 

PiM  +  tijjk  *  *jTfj*<yrf+ 3  Vi/  6  ,*,*  e  *>y 

Pij.*  +  dij,k  <  XijTij,kV;iomyJ 6  V/y  6  6  Pij 

( vj,a )2  -  (Pi,a)~  +  2 (RPpij'C  +  XPqipa  +  Pij^  +  Xjj  qij 'b  +  Rjj  Pij,c  +  Xjj  qij,c) 
>-Mij(\-Xij)  Vij  £  NT 

( Vj,b  f  -  (vub)2  +  2(R^ pij  b  +  xjfqij'b  +  R^j  Pi j,c  +  x’i-qij'C  +  R-j'p^a  +  X^ qijA ) 

>  —  Xjj)  Vij  £  NT 


(2a) 

(2b) 

(3a) 

(3b) 

(4) 

(4a) 

(5a) 

(5b) 
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(5c) 


(vj,cY  -  (Vi,cy  +  2 {RCijPij,c  +  xu<lij,c  +  RjJ Pi j, a  +  X]JqiM  +  R--Pij,b  +  Ku<lij,b) 


rcb , 


>  -M,j(\  —  Xij)  Vij  £  NT 

(Vj,a)2  -  (Vi,a)2  +  2 {R^Pij'd  +  X^a<Jij,a  +  Rjj  Pij,b  +  Xij  +  Ejf  Pij,c  +  X^j  <Jij,c) 
<  Mij{\  -  Xij)  V ij  g  NT 

( vj,b)2  ~  (Vi,by  +  2 {R^Pij'b  +  xffqijj,  +  Rjj  pij,c  +  Xjj  qg ;C  +  R^pg,a  +  Kg' <lij,a) 


<  Mg{\  -Xij)  Vij  £  NT 

( Vj,c )2  _  (vi,c)2  +  2 {RgPij,c  +  XfjC qijtC  +  Rjj  Pi ij>a  +  Xjj  qjj^a  +  RjjjPij,b  +  Xdg (pj,b) 


<  Mjj{  1  -  Xij) 

(0.95 v;wmf  <  vfk  <  (i.05v;iom)2 
0  <  qplM  <  gpui  k 

gquk  <  gq\[k 
-  (R/i.k  <  gq“k 


Vi  j  £  NT 
Vi  g  N,k  e  Vj 
Vi  G  N,  k  G  Pi 
Vi  G  N,  k  G  Pi 
Vi  G  N,  k  G  Pi 


lpi,k  <  dpkk 
lqi,k  —  dqlM 


Vi  &N,k  G  Pj 
Vi  £  N,k  £  Pi 


(5d) 

(5e) 

(5f) 

(6) 

(7a) 

(7b) 

(7c) 

(8a) 

(8b) 


Topology  constraints 

y  xij<\N{o\-\ 

ijeS(C) 

Xij  +  Xjj  <  1 
Xij  +  Xjj  >  1 


VC  (9) 

V(iJ)eS  (10a) 

V(i,j)£&\S  (10b) 


3.1.2  Model  Discussion 

The  objective  function  (1)  minimizes  the  cost  associated  with  all  real  and  reactive  load 
shed  due  to  the  loss  of  generation  capacity  at  bus  i.  The  cost  parameters,  Api,k  and  Aqlk 
(nominally  set  to  1)  scale  the  total  value  of  Kw  and  Kilovolt-ampere  reactive  (Kvar)  shed 
by  the  model. 

Constraints  (2a)  and  (2b)  maintain  balance  of  flow  throughout  the  network,  node  by  node, 
between  generation,  load,  and  load  shedding,  for  real  and  reactive  power,  respectively. 
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Constraints  (3a)  and  (3b)  ensure  that  real  and  reactive  power  flows  for  each  phase  are  within 
stated  tolerances.  A  large  difference  in  flow  on  a  single  line  between  two  phases  isn’t 
desirable  from  an  EE  perspective.  The  fiij  terms  are  tolerance  parameters:  the  higher  the 
Pij  value,  the  wider  the  tolerance  between  phases  that  is  allowed. 

Constraint  (4)  ensures  the  thermal  limits  of  line  (z,  j)  are  not  exceeded.  These  limits  are 
enforced  only  if  a  line  is  being  used  (i.e.,  Xjj  =  1).  Constraint  (4a)  is  a  linear  version  to 
allow  each  x-,j  to  control  whether  flow  is  permitted  on  arc  (z,  j). 

Constraints  (5a)  through  (5f)  account  for  flow  losses  across  all  possible  phase  combinations 
on  line  (z  j)  due  to  resistance  and  reactance  values  for  that  line.  Again,  these  constraints  are 
enforced  only  if  a  line  is  being  used  (i.e.,  x,y  =  1). 

Constraint  (6)  maintains  voltage  at  node  z  within  95%  to  105%  of  the  nominal  voltage  for 
that  node,  as  mapped  by  its  location  from  the  nearest  upstream  transformer. 

Constraints  (7a),  (7b),  and  (7c)  ensure  nodes  only  generate  real  and/or  reactive  power  if 
capable,  as  appropriate. 

Constraints  (8a)  and  (8b)  only  allow  load  to  be  shed  at  nodes  with  real  or  reactive  load. 
This  prevents  the  formulation  from  spuriously  shedding  load  at  nodes  without  load,  simply 
to  satisfy  constraint  mathematics. 

Constraint  (9)  precludes  cycles  in  the  network  by  ensuring  the  number  of  active  edges  in 
a  cycle  is  one  less  than  the  number  of  edges  in  that  cycle.  This  constraint  is  implemented 
differently  than  the  others.  Initially,  we  solve  the  ACPF  model  without  the  cycle  breaking 
constraint  active.  Once  a  solution  is  obtained  we  analyze  it  to  identify  cycles  within  that 
solution.  For  each  cycle  in  the  incumbent  solution,  a  copy  of  constraint  (9)  is  added  and 
then  the  ACPF  model  is  solved  again. 

Constraint  (10a)  ensures  there  is  flow  in  at  most  one  direction  (namely,  on  one  arc)  for  a 
given  edge,  and  constraint  (10b)  ensures  that  edges  without  a  switch  are  active  in  at  least 
one  direction.  If  the  switch  on  edge  (z,  j)  is  open,  then  both  x,j  and  Xjl  will  be  zero. 


15 


3.2  Model  Implementation 

The  model  is  built  in  the  Python  programming  language,  an  interpreted,  high-level,  general 
purpose  programming  language  (Python  Software  Foundation,  2017).  The  modular  and 
flexible  nature  of  Python  enabled  us  to  integrate  all  computational  tasks  within  a  single  pro¬ 
gram.  The  data  read  and  processing  steps  use  several  standard  Python  packages,  including 
j  son,  numpy,  and  math  packages,  as  well  as  networkx  (Hagberg  et  al.,  2008)  for  analysis. 
The  mathematical  program  is  implemented  in  Pyomo,  a  “Python-based,  open-source  op¬ 
timization  modeling  language  with  a  diverse  set  of  optimization  capabilities”  (Hart  et  al., 
2012,  2011). 

Our  test  networks  and  the  real  life  network  are  all  in  JavaScript  Object  Notation  (JSON) 
format  (Ecma  International,  2017),  with  three  dictionaries  in  the  data  file,  one  for  the  nodes, 
one  for  the  edges  in  the  network,  and  another  for  more  detailed  information  on  the  nodes. 
JSON,  consisting  of  name  and  value  pairs  for  its  data  storage  is  a  natural  fit  for  the  Python 
(and  Pyomo)  objects  we  use  to  construct  the  model. 

After  constructing  all  of  the  Python  storage  objects  necessary  for  the  model  we  have  two 
primary  loops  to  prepare  the  data,  one  for  the  nodes  and  one  for  the  edges.  During  these 
loops  we  read  in  all  of  the  requisite  details  for  the  distribution  system.  First,  we  read  in  node 
details  including  attributes  such  as  the  type  of  node,  e.g.,  substation,  shunt  capacitor,  or 
load,  a  unique  node  identification  number,  and  specific  data  items  such  as  dp^k,  dq^k,  g puk, 
vl. , ,  and  more.  Next,  we  read  in  data  for  the  lines  (i.e.,  edges)  in  the  system.  We  capture  the 
Rjj^k  and  for  lines,  which  nodes  in  the  network  a  given  line  connects,  and  the  number 
of  phases  in  a  line.  For  lines  that  are  transformers  we  capture  the  voltage  difference  across 
the  transformer  and  map  these  voltages  through  the  network  as  V[wm. 

We  store  all  data  in  Python  dictionaries.  Dictionaries  are  an  ideal  storage  medium  for  the 
network  data  because  their  structure  uses  key-value  pairs  and  in  this  case  we  use  the  unique 
node  identification  number  to  index  node  data,  and  the  (tail,  head)  node  identification  pairs 
for  an  edge  to  index  edge  data.  Furthermore,  Python  dictionaries  naturally  integrate  with  the 
Pyomo  optimization  package  for  forming  Pyomo  objects.  For  example,  once  a  dictionary 
is  created  with  Rjj,k,  it  can  be  fed  directly  to  the  Pyomo  resistance  object. 

Once  the  data  is  read  in  and  prepared,  we  build  the  ACPF  model  in  Pyomo  and  solve  it  with 
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the  CPLEX  solver  (IBM,  2017). 


3.3  Model  Testing 

We  begin  with  a  basic  test  circuit  that  has  no  cycles  or  switches,  a  minor  design  simplification 
to  facilitate  incremental  model  construction.  A  depiction  of  the  basic  test  circuit  appears  in 
Figure  3.1,  with  associated  parameters  in  Table  3.1  and  Table  3.2. 


Substation  (3,0) 

t 


Substation 

Transformer 


(TO)  (3,5) 


Load  (3,139) 


Service 

Transformer 


(3,119) 


(3,139) 


I 

Shunt  Cap  Bank  (3,139) 


The  Basic  Test  Circuit  consists  of  a  substation  with  functionally  infinite  capacity, 
a  transformer,  connective  line  leading  to  a  service  transformer,  and  a  commer¬ 
cial  grade  customer  with  moderate  loading  and  a  shunt  capacitor  bank.  When 
operated,  no  load  should  be  shed  in  this  test  network. 

Figure  3.1.  Basic  Test  Circuit 


3.4  Model  Performance  on  Basic  Test  Circuit 

The  ACPF  model,  running  on  an  early  2015  MacBook  Pro  with  a  2.7  GHz  i5  processor, 
solves  the  basic  test  circuit  in  less  than  a  second.  The  objective  function  (1)  returns  a  result 
of  0.00  corresponding  to  zero  load  shed,  as  anticipated.  We  observe  pij^  flows  from  the 
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Table  3.1.  Basic  Test  Circuit  Node  Data 


Node 

Type 

9Pi,k 

cKk 

dpi,k 

dqi.k 

1 ynom 

(3,0) 

Substation 

OO 

OO 

0 

0 

69 

(3,5) 

Connection 

0 

0 

0 

0 

7.62 

(3,119) 

Connection 

0 

0 

0 

0 

7.62 

(3,139) 

Load, 

Shunt 

capacitor 

0 

400 

57.55  (a) 
55.48  (b) 
60.65  (c) 

29.29  (a) 
27.51  (b) 
30.78  (c) 

0.277 

Load  for  each  phase  (a-c)  indicated  explicitly 


Table  3.2.  Basic  Test  Circuit  Edge  Data 


Edge 

Type 

Phases 

Tij,k 

Rq,k 

Xii,k 

(3,0), (3, 5) 

Transformer 

a ,  b,  c 

120.773 

0.0217 

0.3548 

(3, 5), (3, 119) 

Line 

a,  b,  c 

498.75 

0.1334  (a) 
0.1354(b) 
0.1343  (c) 

0.3124  (a) 
0.3053  (b) 
0.3089  (c) 

(3,1 19), (3, 139) 

Line 

a,  b,  c 

32.808 

0.0009 

0.0053 

Rjjj,  and  Xij.k  for  each  phase  listed  explicitly  when  variations  exist. 


substation  to  the  load  of  57.55,  55.48,  and  60.65  Kw  for  phases  a,b,  and  c,  respectively, 
values  matching  the  demands  in  node  (3,139). 

Due  to  the  presence  of  a  shunt  capacitor  bank  on  node  (3,139)  all  reactive  demand  is 
satisfied  locally  at  the  commercial  customer,  resulting  in  no  qij^  flow.  Iterative  testing, 
namely  removal  of  the  shunt  capacitor  bank,  causes  qq^  flow  (in  Kvar)  satisfying  all  dqi  k. 
This  behavior  confirms  proper  ACPF  model  operation  regarding  reactive  power  flow  as 
well. 

At  this  point  our  ACPF  model  has  performed  as  expected  on  the  basic  test  circuit,  success¬ 
fully  operating  while  adhering  to  power  flow  constraints  (2a)  through  (8b). 

3.5  Model  Testing  With  Cycles 

Once  satisfied  with  operation  on  the  basic  test  circuit  we  added  switching  controls  and  the 
cycle  breaking  constraint  (9)  to  complete  implementation  of  the  formulation  in  our  ACPF 
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model.  At  this  point  we  build  a  basic  test  circuit  with  a  cycle  and  switches  for  testing,  as 
depicted  in  Figure  3.2,  with  associated  data  in  Table  3.3  and  Table  3.4. 


Switch  2,4 


The  Basic  Test  Circuit  With  Cycle  consists  of  a  substation  with  functionally  infinite 
capacity  and  a  moderate  load,  each  with  their  own  transformer.  Between  them  is 
a  cycle  of  four  nodes,  with  two  of  the  four  edges  containing  switches.  Unlike  the 
Basic  Test  Circuit,  there  is  no  shunt  capacitor  at  the  load. 

Figure  3.2.  Basic  Test  Circuit  with  Cycle 


Table  3.3.  Basic  Test  Circuit  with  Cycle  Node  Data 


Node 

Type 

9Plk 

dpi,k 

d(Ji,k 

ynom 

Sub 

Substation 

oo 

oo 

0 

0 

100 

1 

Connection 

0 

0 

0 

0 

10 

2 

Connection 

0 

0 

0 

0 

10 

3 

Connection 

0 

0 

0 

0 

10 

4 

Connection 

0 

0 

0 

0 

10 

Load 

Load 

0 

0 

30.0 

20.0 

1 

Loading  for  each  phase  a  -  c  is  identical. 
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Table  3.4.  Basic  Test  Circuit  with  Cycle  Edge  Data 


Edge 

Type 

Phases 

Tij,k 

Rij,k 

*ij,k 

Sub,l 

Transformer 

a,  b,  c 

100 

0.0217 

0.3548 

1,2 

Line 

a,  b,  c 

100 

0.1334 

0.3124 

1,3 

Line 

a,  b,  c 

100 

0.1334 

0.3124 

2,4 

Line 

a,  b,  c 

100 

0.1334 

0.3053 

3,4 

Line 

a,  b,  c 

100 

0.1343 

0.3089 

4, Load 

Transformer 

a,  b,  c 

100 

0.0009 

0.0053 

Rjj.k  and  X/jk  for  each  phase  listed  explicitly  when  variations  exist. 


The  Basic  Test  Circuit  with  Cycle  contains  realistic  R/j  and  Xg  values,  drawn  from  real 
networks.  Possessing  a  single  load  with  its  own  service  transformer  and  a  substation  with 
functionally  infinite  gpt,k  and  gqt.k,  the  basic  test  circuit  with  cycle  is  electrically  sound  and 
designed  to  shed  no  load  when  solved. 


3.6  Model  Performance  on  Basic  Test  Circuit  with  Cycle 

Our  model  solves  the  basic  test  circuit  with  cycle  initially,  not  enforcing  the  cycle  breaking 
constraint,  resulting  in  the  behavior  summarized  in  Table  3.5. 


Table  3.5.  Basic  Test  Circuit  with  Cycle  Results  before  Cycle  Breaking 
Constraint 


Edge 

Type 

Switch  /  Position 

xij 

PijX 

Qij-k 

Sub,l 

Transformer 

Lalse,  N/A 

1 

30.0 

20.0 

1,2 

Line 

Lalse,  N/A 

1 

15.0 

10.0 

1,3 

Line 

True,  Closed 

1 

15.0 

10.0 

2,4 

Line 

True,  Closed 

1 

15.0 

10.0 

3,4 

Line 

Lalse,  N/A 

1 

15.0 

10.0 

4, Load 

Transformer 

Lalse,  N/A 

1 

30.0 

20.0 

Pij,k  and  qijtk  flows  the  same  for  all  phases. 
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At  this  point,  we  check  the  incumbent  solution  for  cycles  (using  the  networkx  package) 
and  identify  the  cycle  consisting  of  nodes  1,  2,  3,  and  4.  The  cycle  breaking  constraint  (9) 
is  enforced  and  the  ACPF  model  runs  again,  this  time  producing  a  solution  with  no  cycles 
and  zero  load  shed,  as  summarized  in  Table  3.6.  Cycle  breaking  is  restricted  to  opening  a 
switch  on  lines  that  have  switches. 


Table  3.6.  Basic  Test  Circuit  with  Cycle  Results  after  Cycle  Breaking  Con¬ 
straint 


Edge 

Type 

Switch  /  Position 

xij 

Pij,k 

Qij,k 

Sub,l 

Transformer 

False,  N/A 

1 

30.0 

20.0 

1.2 

Line 

False,  N/A 

1 

0.0 

0.0 

1,3 

Line 

True,  Closed 

1 

30.0 

20.0 

2,4 

Line 

True,  Open 

0 

0.0 

0.0 

3,4 

Line 

False,  N/A 

1 

30.0 

20.0 

4, Load 

Transformer 

False,  N/A 

1 

30.0 

20.0 

Pij,k  and  qtj'k  flows  the  same  for  all  phases. 


As  can  be  observed,  the  cycle  breaking  constraint  has  forced  the  switch  on  line  2,4  to  open. 
Flow  redirected  to  meet  demand  while  avoiding  a  cycle  being  present  in  the  network. 

3.7  Model  Validation 

Having  tested  our  ACPF  model’s  basic  functions  we  test  it  against  an  IEEE  RTS  network  for 
validation.  The  IEEE  RTS  is  a  standardized  test  circuit  whose  behavior  is  well  understood 
by  researchers  and  academia.  Furthermore,  phase  changes,  multiple  branching  paths,  and 
uneven  loading  make  it  a  more  robust  test  case  than  the  basic  test  networks  used  thus  far. 
Observing  proper  behavior  of  the  ACPF  model  on  the  RTS  network  serves  as  validation  to 
our  modeling  and  these  efforts  are  described  in  Chapter  4. 
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CHAPTER  4: 
Results 


4.1  IEEE  RTS  13-bus  feeder 

We  choose  the  IEEE  RTS  13  node  test  feeder  because  it  is  relatively  small  and  frequently 
used  to  test  analysis  software  (IEEE  Power  and  Energy  Society,  2017).  Furthermore,  the 
network  itself  mimics  real  world  distribution  and  is  “characterized  by  being  short,  relatively 
highly  loaded,  a  single  voltage  regulator  at  the  substation,  overhead  and  underground  lines, 
shunt  capacitors,  an  in-line  transformer,  and  unbalanced  loading”  (IEEE  Power  and  Energy 
Society,  2017).  A  network  with  these  characteristics  will  test  many  facets  of  our  model 
in  a  more  robust  capacity  than  the  test  networks  in  Chapter  3.  The  network  is  constructed 
without  the  potential  for  cycles  however,  so  our  cycle  breaking  constraint  is  not  exercised 
when  solving  the  IEEE  RTS  13  node  test  feeder.  The  circuit  is  depicted  in  Figure  4.1. 


Figure  4.1.  IEEE  13  Node  Test  Feeder.  The  substation  is  at  node  634 
with  its  transformer  between  node  634  and  633.  A  voltage  regulator  exists 
between  node  632  and  650  but  has  no  impact  on  network  behavior  within 
our  model  because  node  650  is  an  unloaded  bus. 
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Full  electrical  characteristics  for  the  13  node  test  feeder  can  be  found  on  the  IEEE  Power  and 
Energy  Society  website  (IEEE  Power  and  Energy  Society,  2017).  The  basic  characteristics 
and  those  specifically  pertinent  to  our  ACPF  model  are  described  in  Tables  4.1  and  4.2. 
Additionally,  there  is  a  distributed  load  in  the  IEEE  13  node  test  feeder  along  the  line 
between  nodes  632  and  671.  We  capture  this  load  as  an  extra  node  named  “dload,”  existing 
equidistant  between  nodes  632  and  671.  The  distributed  load  node  has  dpkk  values  of  10, 
38,  and  68,  with  dq kk  values  of  17,  66,  and  1 17,  all  for  phases  a,  b,  and  c,  respectively. 


Table  4.1.  IEEE  13  Node  Test  Feeder  Node  Data 


Node 

Type 

gplk 

9<tLk 

dpLk 

dcp.k 

ynom 

634 

Substation, 

Load 

oo 

oo 

160(a) 

120(b,c) 

110(a) 

90(b,c) 

115 

633 

Connection 

0 

0 

0 

0 

4.16 

632 

Connection 

0 

0 

0 

0 

4.16 

650 

Connection 

0 

0 

0 

0 

4.16 

645 

Load 

0 

0 

170(b) 

125(b) 

4.16 

646 

Load 

0 

0 

230(b) 

132(b) 

4.16 

671 

Load 

0 

0 

385(a,b,c) 

220(a,b,c) 

4.16 

684 

Connection 

0 

0 

0 

0 

4.16 

652 

Load 

0 

0 

128(a) 

86(a) 

4.16 

611 

Load, 

Shunt 

capacitor 

0 

100(c) 

170(c) 

80(c) 

4.16 

680 

Connection 

0 

0 

0 

0 

4.16 

692 

Load 

0 

0 

170(c) 

151(c) 

4.16 

675 

Load, 

Shunt 

capacitor 

0 

200 

(a,b,c) 

485(a) 

68(b), 290(c) 

190(a) 

60(b), 2 12(c) 

4.16 

Load  for  each  phase  (a-c)  indicated  explicitly 
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Table  4.2.  IEEE  13  Node  Test  Feeder  Edge  Data 


Edge 

Type 

Phases 

Tij,k 

634,633 

Transformer 

a,  b,  c 

1500 

633,632 

Line 

a,  b,  c 

340 

632,650 

Line 

a,  b,  c 

730 

632,645 

Line 

b,  c 

230 

645,646 

Line 

b,  c 

230 

632,dload 

Line 

a,  b,  c 

730 

dload,671 

Line 

a,  b,  c 

730 

671,684 

Line 

a,  c 

230 

684,611 

Line 

c 

230 

684,652 

Line 

a 

165 

671,680 

Line 

a,  b,  c 

730 

671,692 

Line 

a,  b,  c 

1500 

692,675 

Line 

a,  b,  c 

260 

Rjj  j,  and  Xjj  j.  for  the  network  are  3  x  3  matrices.  We  refer  our  readers  to  the  IEEE  RTS 

reference  materials  for  specific  values. 


The  network  structure  and  parameters  are  built  into  the  JSON  format  our  associates  at 
LANL  are  using  and  which  our  model  is  built  to  read. 


4.2  IEEE  13  Node  Test  Feeder  Performance 

The  IEEE  RTS  13  node  test  feeder  runs  in  a  little  over  a  second  and  returns  an  optimal 
solution  objective  function  value  of  0.00,  a  result  our  associates  at  LANL  confirmed. 
Electrical  behavior  is  captured  in  Table  4.3. 
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Table  4.3.  IEEE  13  Node  Test  Feeder  Edge  Results 


Edge 

Type 

Xij 

Pij,k 

Qij,k 

634,633 

Transformer 

1 

1015  (a) 

381.91  (a) 

919  (b) 

515  (b) 

1132  (c) 

451  (c) 

633,632 

Line 

1 

1015  (a) 

381.91  (a) 

919  (b) 

515  (b) 

1132  (c) 

451  (c) 

632,650 

Line 

1 

0.00 

0.00 

632,645 

Line 

1 

400  (b) 

257  (b) 

645,646 

Line 

1 

230  (b) 

132  (b) 

632,dload 

Line 

1 

1015  (a) 

381.91  (a) 

519(b) 

258  (b) 

1132  (c) 

451  (c) 

dload,671 

Line 

1 

998  (a) 

371.91  (a) 

453  (b) 

220  (b) 

1015  (c) 

383  (c) 

671,684 

Line 

1 

128  (a) 

86  (a) 

170  (c) 

0.00  (c) 

684,611 

Line 

1 

170  (c) 

0.00  (c) 

684,652 

Line 

1 

128  (a) 

86  (a) 

671,680 

Line 

1 

0.00 

0.00 

671,692 

Line 

1 

485  (a) 

115.11  (a) 

68  (b) 

21.33  (b) 

460  (c) 

197  (c) 

692,675 

Line 

1 

485  (a) 

65.91  (a) 

68  (b) 

0.00  (b) 

290  (c) 

12(c) 

Our  ACPF  model  demonstrates  full  functionality  at  this  point,  managing  cycles,  optimal 
switching,  and  returning  expected  operational  results  on  a  standardized  test  circuit. 
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Additionally,  running  the  ACPF  model  against  variations  in  our  test  networks  (e.g.,  removing 
components,  constraining  electrical  parameters,  and  varying  loading  amounts)  produced 
electrical  results  we  would  expect  for  these  variations. 


4.3  Expansion  of  IEEE  13  Node  Test  Feeder 

Our  model  is  ultimately  about  performing  analysis  on  distribution  systems  under  the  stress 
of  disruption.  With  that  in  mind,  we  build  an  expanded  version  of  the  IEEE  13  node  test 
feeder,  essentially  cloning  the  network,  with  “a”  and  “b”  versions  of  each  node  and  line. 
The  two  versions  connect  via  a  line  from  node  671a  to  671b.  This  3-phase  line  contains  a 
switch  which  is  open  under  normal  operating  conditions.  This  expanded  IEEE  13  node  test 
feeder  is  depicted  in  Figure  4.2 


The  duplicated  13  node  test  feeder  with  “a”  and  “b”  sides,  connected  by  a  line 
between  nodes  671a  and  671b.  Electrical  parameters  for  both  sides  are  the  same 
as  the  original  IEEE  13  node  test  feeder,  described  in  Tables  4.1  and  4.2.  The 
cross  connection  line  is  3-phase  with  Rjj  and  Xjj  values  copied  from  line  692-675. 
The  Tij  value  is  very  low  at  15  amps. 

Figure  4.2.  Duplicated  13  Node  Test  Feeder 


When  given  to  our  ACPF  model  the  duplicated  IEEE  13  node  test  feeder  sheds  no  load, 
an  expected  result  since  the  single  version  met  all  load  demands.  In  essence,  we  have 
parallel  instances  of  the  previous  test  network  and  each  one  saw  commensurate  electrical 
behavior,  analogous  to  the  results  in  Table  4.3  for  both  the  “a”  and  “b”  sides.  Our  next 
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step  is  to  consider  model  performance  and  distribution  system  behavior  in  the  presence  of 
a  disruption. 


Our  chosen  disruption  is  failure  of  the  “b”  side  substation,  essentially  removing  it  from  the 
network.  Running  the  model  with  this  disruption  results  in  all  “a”  side  load  satisfied  and  the 
cross  connection  switch  to  close.  Including  the  cross-connection,  only  a  handful  of  lines 
on  the  “b”  side  have  flow.  Their  values  are  captured  in  Table  4.4. 


Table  4.4.  Duplicated  IEEE  13  Node  Test  Feeder  Edge  Results 


Edge 

Type 

xij 

Pij,k 

Qij,k 

67  la, 67  lb 

Line 

1 

58.07 

58.07 

671b,dload_b 

Line 

1 

9.52  (a) 

38.82  (b) 

41.40(c) 

5.29  (a) 

13.90  (b) 

18.24  (c) 

Pij,k  and  %,k  A 

lows  the  same  for  all  phases  unless  exp] 

icitly  listed. 

Due  to  relatively  low  thermal  limit  we  imposed  on  the  671a  -  671b  connection,  most  of  the 
“b”  side  goes  without  power  and  the  model  returns  an  objective  function  value  of  4689.38. 
We  are  not  surprised  by  this  amount  of  load  shed  given  the  parameters  we  assigned  the 
cross-connection. 


4.4  Model  Performance  Overview 

At  this  point  the  model  has  demonstrates  the  capability  to  analyze  standardized  test  circuits 
in  the  presence  of  network  disruptions.  Electrical  parameters  behave  as  expected  and 
operation  of  any  distribution  system  can  be  optimized  by  the  ACPF  model  in  an  effort  to 
minimize  load  shed.  We  provide  recommendations  for  follow-on  work  and  additional  model 
applications  in  Chapter  5. 
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CHAPTER  5: 

Conclusions  and  Follow-on  Work 


5.1  Conclusions 

Working  within  a  Python  and  Pyomo  environment  with  CPLEX  as  a  solver,  our  model 
demonstrates  the  ability  to  find  an  optimal  network  configuration  and  associated  power  flow 
for  a  3-phase  electrical  distribution  system.  Our  ACPF  model  can  handle  phase  changes, 
varying  and  imbalanced  loads,  real  and  reactive  demand,  switch  operation,  and  ensures  the 
network  is  configured  as  a  spanning  tree  to  meet  standard  distribution  system  behavior. 

The  novelty  in  our  work  is  due  to  the  focus  on  distribution  systems.  The  model  receives 
input  from  a  JSON  data  structure  used  by  the  individuals  at  LANL  working  on  infrastructure 
analysis.  This  allows  for  downstream  use.  Furthermore,  the  free  and  open-source  nature  of 
the  Python  programming  language  makes  our  work  easily  accessible  for  additional  research. 

Alternating  current  is  inherently  non-linear.  That  said,  the  mathematics  of  the  model  are 
linearized,  with  the  exception  of  constraint  (4)  which  is  non-linear.  Testing  and  experi¬ 
mentation  revealed  that  the  non-linear  structure  of  constraint  (4)  allowed  for  nonsensical 
solutions  when  implemented  independently.  As  such,  we  identified  the  need  to  build  and 
include  constraint  (4a).  Doing  so  achieves  correct  results  from  the  test  networks  and  the 
IEEE  RTS  13  node  test  feeder. 


5.2  Follow-on  Work 

Our  associates  at  LANL  have  a  detailed  and  highly  comprehensive  dataset  from  a  real 
northeastern  utility  provider.  Our  test  networks  and  IEEE  RTS  13  node  test  feeder  are  built 
in  the  same  JSON  format  as  the  real-world  data,  facilitating  model  application  to  the  real 
world  data  after  some  preliminary  data  cleaning. 

The  real-world  dataset  contains  52  substations,  73,304  loads,  and  several  hundred  thousand 
nodes  and  arcs.  There  is  interest  in  systematically  triggering  failure  of  substations  and 
capturing  network  behavior  after  the  ACPF  model  runs.  This  analysis  would  yield  data  for 
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subsequent  vulnerability  analysis  as  test  and  training  sets.  Furthermore,  it  would  inform 
us  regarding  distribution  system  behavior  in  the  presence  of  substation  failure,  a  topic  of 
interest  to  various  government  agencies. 

To  facilitate  future  work  investigating  the  behavior  of  distribution  systems  in  the  presence  of 
adverse  events,  the  model  should  be  augmented  to  more  easily  disable  part  of  the  distribution 
system  being  analyzed.  At  present,  the  data  itself  needs  to  be  modified  to  break  a  line  or 
disable  a  substation.  An  update  with  the  ability  to  pass  in  a  component  identification  number 
as  a  function  argument  for  disruption  purposes  would  greatly  simplify  this  process. 

Lastly,  our  ACPF  model  is  built  as  a  single  Python  function,  allowing  for  easy  incorporation 
into  AD  and  DAD  analysis.  Analysis  of  this  nature  would  yield  insight  to  the  system 
vulnerabilities  of  any  distribution  network  and  the  corresponding  hardening  efforts  of  the 
greatest  benefit. 
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